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We investigate the properties of an atom under the influence of a synthetic three-dimensional spin- 
orbit coupling (Weyl coupling) in the presence of a harmonic trap. The conservation of total angular 
momentum provides a numerically efficient scheme for finding the spectrum and eigenfunctions of 
the system. We show that at large spin-orbit coupling the system undergoes dimensional reduction 
from three to one dimension at low energies, and the spectrum is approximately Landau level-like. 
At high energies, the spectrum is approximately given by the three-dimensional isotropic harmonic 
oscillator. We explore the properties of the ground state in both position and momentum space. 
We find the ground state has spin textures with oscillations set by the spin-orbit length scale. 



I. INTRODUCTION 

The recent experimental success in simulating spin- 
orbit coupling|TH5] in a cold atom context has produced 
great interest in the field. Although the experimental 
setup produced only an Abelian spin-orbit coupling, it is 
possible to produce more complicated non-Abelian syn- 
thetic fields, such as any combination of Rashba and lin- 
ear Dresselhaus couplings. [SHIT] The isotropic limit, such 
as Rashba or linear Dresselhaus coupling, is particularly 
interesting. In this limit, and the absence of an external 
potential, the ground state of the spin-orbit coupling has 
an infinitely degenerate ground state manifold, which is 
given by a ring in momentum space. The ring of min- 
ima gives an energetically free direction for low-energy 
quantum fluctuations. This leads to an effective dimen- 
sional reduction for low energy properties |12IfbI] . such 
as a Landau- like spectrum |12[ HHITfJI , enhanced binding 
energy [L7fE5] . or spontaneous symmetry breaking[26j. 
In the presence of a trap or interactions this degeneracy 
will be broken, up to the two-fold degeneracy guaranteed 
by time-reversal symmetry. However, as we show here, 
effects of the infinite manifold of states will survive even 
to the trapped regime, and the dimensional reduction of 
low energy states will still be visible. 

A recent proposal [161 |2"T1 |2"5] allows for the study 
of spin-orbit coupling that has no solid state counter- 
part: three-dimensional spin-orbit coupling, or Weyl 
coupling[29j, of the form: vp ■ cr, where p is momentum, 
and cr is the spin operator. The conceptual picture in this 
system is similar to that of the Rashba case. However, 
instead of a ring of ground states we now have a spheri- 
cal manifold of ground states. This implies there are two 
energetically free directions for low energy quantum fluc- 
tuations, so the dimensional reduction is now from D = 3 
to D = 1. 

In this paper we analyze the motion of a particle mov- 
ing in an isotropic harmonic trap, subject to Weyl cou- 
pling. We show that while Weyl coupling is the natural 
extension of Rashba coupling to incorporate the third 



spatial dimension, this extension actually increases the 
overall symmetry of the system and simplifies its treat- 
ment. We describe a numerical scheme for finding the 
eigenvalues exactly, and we analyze the eigenfunctions. 
Even the ground state of this system exhibits significant 
spin and orbital textures. In the limit of large coupling 
we find Landau like behavior of the spectrum, in that the 
low-energy spectrum tends to that of a one-dimensional 
harmonic oscillator with a large number of nearly degen- 
erate levels. However, at sufficiently large energies, the 
spin-orbit coupling becomes a perturbation, and the high 
energy spectrum is well approximated by the 3D isotropic 
harmonic oscillator. 

The paper is organized as follows. We first review the 
model Weyl coupling Hamiltonian. We perform a partial 
wave expansion of the system in the presence of a trap, 
and find a numerically efficient method for diagonaliz- 
ing the system. We calculate the spectrum as a func- 
tion of the spin-orbit coupling strength, and find the low 
energy spectrum becomes Landau-like at large coupling. 
At sufficiently high energies, the spectrum crosses over to 
that of the three-dimensional isotropic harmonic oscilla- 
tor. We derive a system of coupled differential equations 
in the radial coordinate only, and use it to explain the 
crossover between the Landau spectrum and isotropic os- 
cillator spectrum. Finally, we show the ground state is 
characterized by persistent orbital and spin currents. 



II. MODEL 

Our Hamiltonian is given by 
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n r 

H=^ + vp-a + -, (1) 

where r and p are respectively the position and momen- 
tum operators and cr — (oi, oi, 03) is the vector of con- 
ventional Pauli matrices. This Hamiltonian can be im- 
plemented in cold atomic systems using two-photon tran- 
sitions as described in |27| : here v is a coupling constant 
subject to external control. The two dimensions of the 
spin operator cr correspond to two internal states of the 
atom, which would ordinarily be two different hyperfine 
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states. We have chosen a system of units in which the 
reduced Planck constant h, the mass M of the particle, 
and the harmonic oscillator frequency lj are all equal to 
1. 

Our Hamiltonian commutes with an angular momen- 



tum J 



where L = r x p is the orbital angular 



momentum of the atomic center of mass. Note that this 
J need not be identical to the usual angular momentum 
that is constructed from the sum of the orbital angular 
momentum of the atom's center of mass and the angular 
momenta of its electrons and nuclei. However, it is a con- 
served quantity whose three spatial components satisfy 
the usual commutation relations of an angular momen- 
tum operator, and in this sense we can treat J as the sum 
of an integer-valued orbital angular momentum L and a 
spin s = h. The er • p term is a scalar under the rota- 
tions induced by J, as is the kinetic and trapping term, 
so our Hamiltonian commutes with J. We note in pass- 
ing that Rashba coupling is obtained from removing the 
<r 3 p 3 term from Weyl coupling. Thus Rashba coupling 
is a linear combination of scalar and rank-2 tensor op- 
erators under rotations, whereas Weyl coupling has the 
simpler scalar form and is much simpler to solve with a 
standard partial-wave decomposition. 

In the absence of a trap, the spectrum of the Weyl 
coupling can be found exactly. There are two bands cor- 
responding to spin aligned, and anti-aligned with mo- 
mentum. The spectrum is given by 



E(P) 



±v\p\ 



(2) 



The low energy band has a minimum defined on the 
sphere |p| = v. Near the minimum of this sphere, the dis- 
persion is parabolic only along the radial direction, and 
is constant along the polar and azimuthal axes. This sug- 
gests that for low energy properties of the system, quan- 
tum fluctuations will be energetic only along the radial 
direction. We expect that in the presence of a spher- 
ically symmetric trap, the low energy spectrum of the 
system will undergo a dimensional reduction from D = 3 
to D = 1 for a sufficiently large spin-orbit coupling. 

The conserved total angular momentum J allows us to 
use a basis of well defined angular momentum \j, m, A), 
where 



A. Number Basis 

The Weyl coupling Hamiltonian can be numerically di- 
agonalized by performing a partial wave decomposition 
into states of the three-dimensional isotropic harmonic 
oscillator, with angular momentum state \j, m,A). We 
define the basis \n,j, m,A) = \n)\j,m,X) in Appendix 
|VI| as a state with n radial quantum nodes, and an an- 
gular momentum eigenstate given by p|. When v = 
these states have energy E = 2n + I + . [3U] It is the 
convenient to express the spin-orbit coupling in terms of 
creation(annihilation) operators a^(a) as 



p • er = —= ycr ■ a — er • aj 

v2 



V2 



(A + - X 



(4) 
(5) 



where A + = <r • a' and A~ = er • a are two rank-0 tensors 
with respect to rotations generated by J. These operators 
satisfy the commutation relation 1 {A + , A~ } = a^ • a + 
| = E/huj, which allows us to express our spin-orbit 
coupled Hamiltonian as 



H=\{A+,A-} 



v^(A+-A- 



(6) 



In Appendix | VI| we show that the matrix elements of the 
operators A + and A~ are given by 



A \n,j,m,-} = \/2(n + j + l)\n,j,m,+) (7) 
A~\nJ,m,+) = V2n\n - m, -) 



and 



A + \n,j,m,+) = y/2{n + j + l)\n,j,m,-) (8) 

A + \n,j,m,-) = y/2(n + l)\n + l,j,m,+). 

Since the operator A + contains a combination of cre- 
ation operators, it is clear that it will raise the energy of 
a state by one unit. In the radial basis there are two ways 
to raise the energy by one unit. The angular quantum 
number can be increased by one with the radial quantum 
number held constant: AZ = +1, An = 0, or the radial 
number can be increases by one and the angular quan- 
tum number lowered by one. Repeated applications of 
the raising operator alternate between state with I = j + s 
and I = j — s, while increasing the energy by one unit 
each time. 



| j, m, A) = e 



j-Xm+xx 

j+Xm+xx 
2(j+xx) 



\3 + 



(3) 



and A = ±1 corresponds to J and er aligned (j = I + s) 
or anti-aligned (j = I — s). These states have eigen- 
values J 2 \ji Tn, A) = j(j + l)\j,m, A) and J z \j, m, A) = 
m\j, m, A). They are complete, which allows us to project 
our Hamiltonian into subspaces of fixed j, m. 



B. Numerical Diagonalization 

The spin-orbit coupled Hamiltonian can therefore be 
numerically diagonalized efficiently by first projecting the 
Hamiltonian into sectors of good \ j, m) 



H = J2Y1 H im\3m)(jm\ 



(9) 
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Figure 1. Numerical calculation of the spectrum as a func- 
tion of v. At v = 0, all energy levels are plotted have en- 
ergy E < 15 + §. As the spin-orbit strength is increased, 
the levels split off into groups with increasing radial quantum 
number n. Each landau level has energy of approximately 

E n jm = (n+ |) + 2l) l V- As a visual guide, each level 

is shifted by « 2 /2. Two regimes are clearly identifiable, corre- 
sponding to a 3D harmonic oscillator with slight level mixing, 
and the ID Landau-level problem. As discussed in the text, 
the crossover between these two regimes is given by E ~ v 2 /2. 
This crossover is shown by the black line. 



where the Hamiltonian Hj m is defined by the matrix el- 
ements 

(n'X'\H jm \n\) = (2n + j + l)*„, n ^ A ,v (10) 
+ iv (Vn + I8 n+I . n ,5 y+ 6 X - 

+ iv^/n+j + ld n ,n' (8\+S\'- - S\'+S\-) ■ 

This matrix is tridiagonal, and can be efficiently di- 
agonalized through 0(n) operations. Figure 1 shows 
the spectrum as a function of the spin-orbit parame- 
ter v. At v = 0, we have included all states that have 
N < 10, where N — 2n + I is the total quanta of the 
3D isotropic oscillator. Two regions in the spectrum are 
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identifiable. For N 3> \ ^ ne spectrum is approximately 
that of the three-dimensional harmonic oscillator with 



by 



^ . For TV <C ^ the spectrum is given 



E 



0+R 

2v 2 



v 
~2 



(11) 



as will be shown in Sec II C 3 To lowest order in 1/v 



this has the form of a one dimensional harmonic oscil- 
lator in the radial mode. We will later see that these 
states are well localized in the momentum space poten- 
tial near |p| = v. This suggests that the energetically 
free excitations along the polar and azimuthal directions 
are not important for low energy states. The system 
therefore undergoes a dimensional reduction from D = 3 
to D = 1. The one-dimensional structure is reminiscent 
of the Landau levels. Since the problem is spherically 



Figure 2. Momentum space eigenf unctions /^(p) for low 
energy states at large spin-orbit coupling. The red, blue 
and black curves correspond to states with n = 1,2,3 re- 
spectively. The eigenfunctions are well approximated by the 
one dimensional harmonic oscillator wavefunction centered 
around p — v. The number of nodes corresponds to the ra- 
dial quantum number. States with higher j have the same 
form, up to a possible sign, provided the number of radial 



quanta are smaller than ~ v /2. 
have opposite sign. 



The eigenfunctions / nj (p) 



symmetric, mixing of the Landau levels appears only at 
higher order in the inverse spin-orbit coupling parameter. 
Note that while these levels have been seen in previous 
work, |121 [T3HHI| , the crossover to the three-dimensional 
spectrum was missed. 



C. The Schrodinger Equation as a system of 
coupled differential equations 

1. Momentum Space 

The eigenstates of the Hamiltonians Hj m are states of 
good total angular momentum. In general they can be 
expressed as 



,J,m) 



m, - 



\i>n r ) li' 



(12) 



where ) is an eigenstate of Hj m with n r radial modes. 
This form suggests that each Hamiltonian Hj m has a 
corresponding set differential equation in only the radial 
degrees of freedom. It is more natural to work in momen- 
tum space, where p is a dynamical variable, instead of an 
operator. The presence of the harmonic trap allows us 
to treat the position operator as a derivative, r = iV p . 
The corresponding "dual" Schrodinger equation is 
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vcr ■ p 



^(p) = £ty(p). (13) 



The radial eigenfunctions in momentum space have the 
corresponding form 

^nim(p) = fnj(p)xJ m (P) + fnjip)Xj m (P), ( U ) 
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where p 
spinors x 



:i>» 



|p|, P = P/P; fnjiP) = (P|^n r )i and thc 
— (p|j m= ' = ) have total momentum j, with 



a J3 projection m. As shown in the appendix, the action 
of the operator <x • p operating on the spinors xf m (p) 



is to interchange the spinors and multiply the result by 

p, i.e., a ■ px^ m (p) = PxJ m (p)- The form (Pf S ives a 
consistent set of two coupled differential equations in the 
independent variable p. These coupled differential equa- 
tions takes the form 



'2dp* 



(i+j) 

2 p 2 



2 dp 2 2 p 2 



2 
2 



(p) + «pu ni b) = 



raj 1 



M rTj (P) + V P U nj 0) = B M nj (?) 



(15) 
(16) 



where u^Ap) 



pftM 



2. Position space 



The eigenfunctions f„j(p) can be found by solving 
these coupled differential equations. Alternatively, they 
can be constructed using the eigenvectors found from di- 
agonalizing the projected Hamiltonians Hj m . Some ex- 
amples of these functions are shown in Fig(2] for a large 
v. The functions are well approximated by the one- 
dimensional harmonic oscillator wavefunctions centered 
around p = v. 



A radial differential equation can be found in position 
space in a manner analogous to the momentum space 
differential equations. The position space analogue of 
( f 4 1 is given by 



(17) 



where 



r A> 9nj( r ) = ( 



r\^ni) and xf m (r) 



(r\jm\). In the appendix we show that the functions 
9ni ( r ) satisfy the radial differential equation 



U^r 
r \dr 2 

1 / d 2 



r 



dr 2 



(J + I)(J + 1) , 



9nj( r ) 



d 

dr 
d 
dr 



3 + 



9 nj (r) 



E nj9nj{r) 



9nj( r ) = E nj9m( r ) 



(18) 
(19) 



These differential equations can be solved to find the 
eigenfunctions in position space. Alternatively, we can 
the radial eigenfunction from the partial wave expansion 
found in Sec|nU 



3. Landau Levels 



The Landau levels we found numerically, which are 



described by (111, can be understood by considering thc 



asymptotic form of the radial differential equations, ( 15 ) 
and (16), in momentum space. At large momentum, 
p > v, the 1/p 2 term becomes negligible, and the differ- 
ential equations can be decoupled by taking an even/odd 
superposition of (151 and (16 1. The two decoupled dif- 



ferential equations are 

1_ £P_ - 2 
2~df- ' 



± vp 



U nj(P) = Eii 



"J 



(P). (20) 



where = l(it + ±u~). These differential equations 
can be mapped to the one-dimensional harmonic oscilla- 
tor by performing a change of variables p-yp^w. The 

solutions are given by u^{p) = c n H n (p ± t>)e _<P 2 ' , 
where H n (p) is the n-th Hermite polynomial, [31_ and c„ 
is a normalization constant. However, the equation for 
u + (p) is localized around p = —v, outside of the range of 
definition of the radial coordinate. We therefore assume 
these solutions are uZAp) = 0. Transforming, we find 



1 1 

2™n! tt 1 /* 



H n (p-v)e- {p - v ^/ 2 



that u + (p) = —u (p) = 
in the asymptotic limit. 

The form of the harmonic oscillator Hamiltonian sug- 
gests the spectrum is 

1 v 2 

E n =n+~-—. (21) 

To lowest order in v, this result is consistent with the 
large Landau- level like degeneracy found earlier. This de- 
generacy is lifted by the centrifugal barrier, which mixes 
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(a) Ground State 
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(b) Lowest radial numbers 




wavefunction, this energy shift is 



SK ■ = - 



3 + 



3 + 



2 n nU 



dp 



P 



(23) 



Formally, the integral in ( 23 1 is divergent as p — > 0. How- 



ever, we can cure this by introducing a low energy cutoff, 
and subtracting the divergent contribution. The integral 
is then dominated near p = v, and can be well approx- 

imated by ^ dp ^^'^' = ^+OH 



which is valid as long as (H n (p — v)) 2 e ^ is local- 
ized in a region away from p = 0. The lowest order shift 

(j+iY 

j njm — 2 1 ; an d the spectrum in the 



in energies is 



6E„ 



asymptotic limit is 



(3 + 



2v 2 



v 

y 



(24) 



0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 



consistent with the spectrum found using numerical di- 
agonalization in the previous section. 



(c) Lowest j values 




4- Validity of Approximation 



Figure 3. (a) Position space eigenfunctions for low ground 
states at large spin-orbit coupling. The red and blue curves 
correspond to the functions rg^(r) and rg~^(r) respectively. 

(b) The first three position space eigenfunctions rg^(r) for 
j — 1/2 and n — 1, 2, 3 (red, blue and black respectively.) (c) 
The position space eigenfunctions rg„j (r) in the lowest radial 
mode for three lowest values of j = !,§,§, correspond to the 
red, blue and black curves respectively. 



the Landau levels near p = 0. To lowest order in pertur- 
bation theory, the energy shift in the state ipnjm is given 

by 

u 2 



I P 

where I 2 is the orbital angular momentum operator. Us- 



To derive ( 24 1 , we considered our particle in a combi- 



ing ( 14 1 and the asymptotic expression for the radial 



nation of two "spin'-dependent potentials. The first is the 
centrifugal barrier, and the second is the spin-orbit cou- 
pling. The centrifugal barrier is repulsive, and divergent 
at p = 0. The spin-orbit term will form a well centered 
at p = v. For sufficiently large barriers, the minimum of 
the well will be far from the region where the centrifugal 
barrier is finite. This implies that low energy states will 
be well localized in the potential minimum produced by 
the spin-orbit coupling. Since the well is approximately 
harmonic near the minimum, the radial wavefunctions 
will be given by the one-dimensional harmonic oscillator. 
These states are exponentially localized near p = v, and 
will have minimum overlap with the centrifugal potential. 

Wavefunctions with larger radial quantum numbers 
will be increasingly delocalized. The centrifugal barrier 
cannot be neglected when the momentum-space wave- 
function is finite near p = 0. In this limit, the two states 
with j — l±s become mixed with the centrifugal barrier. 
For higher radial quantum numbers, the system is better 
described by two radial harmonic oscillators with angular 
momentum I = j ±s, and the spin-orbit term acts to mix 
the two states. Thus, even when v 3> 1, there will be 
a critical atomic number for which higher energy states 
have average momentum (p) ~ 2v, where the p 2 /2 ki- 
netic term dominates the —vp of the spin-orbit coupling. 
Above this threshold, the effect of the spin-orbit cou- 
pling will become a perturbation, and the spectrum will 
approximate the three-dimensional harmonic oscillator. 
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III. GROUND STATE 

In the previous section we showed the momentum 
space ground state of a trapped particle with Weyl cou- 
pling has the approximate form 



(a) Spin current 



(b) Orbital current 



</>0 A±i(p) = 



1 

71 



Xi 



i(p) -xt + i(i>] 



-(p-vf/2 



7TV4 

(25) 

for v ^> 1. The position space wavefunction 
can be found using the radial differential equa- 
tions. Alternatively, we can directly apply the 
Fourier transform to the momentum-space wavefunc- 
tion, Vo( r ) = / d 3 pe tp ' r ipo(p)- This is most eas- 
ily evaluated by expanding the exponent e lpr — 
^12imi l 3i(P r )Yi m (*) (Xi m (p))* ■ Integration over the 
angular coordinate converts the spinor xf m (p) *° xf m (^)- 
The radial component is then found from the integral 

J °° p 2 ji±i(pr) e ^i/i dp. The exponential factor lo- 
calizes the integrand to a region near p ~ v, in the asmyp- 
totic regime with v 3> 1, this integral can be evaluated by 
extending the lower limit to — oo, and then using the ex- 
plicit form of the spherical Bessel functions, jo(x) = 
and ji(x) — — The position-space radial 

wavefunctions for the ground state are therefore given by 

V2 1 



/o(r) = 



r 3/4 
V2 



2 /o 

(r co$(rv) + v sin(rv )) e~ r 1 



(26) 



V 2 1 2 
/i( r ) = ~^77^> ((1 + t 2 ) sin(rw) — rvcos(rv)) e~ r 0,7) 

with asymptotic corrections of 0(l/v 2 ). The full position 
space wavefunction is approximately given by 

a ±i0) = 4= f/oWxT + i(f)- i/iWxt + i(f) 



72 

in the asymptotic limit. 



(28) 



A. Spin Textures and Currents 

1. Momentum Space 
The azimuthal component of the spinor xti suggests 

2 2 

the Weyl coupling should have a non-zero current in 
the ground state. We first define the current opera- 
tor in momentum space in the standard way. First 
we multiply Schrodinger's equation by ip\ —itp^-^tp = 
r 2 1 

ijp — f- + ip + vp- Njt<Tip\ , and subtract the result 
from it's complex conjugate. The momentum-space con- 
tinuity equation is 

(29) 



where 



[^V p ^ - V T V P V*] = 3 [^VpV] (30) 





(c) Total current 




Figure 4. Spin, orbital and total currents for the ground state 
of a trapped particle with Weyl coupling, with j = i, m = | 
and v — 5. In all figures the plane defined by y — plotted. 
The arrows represent flows of the local normalized current vec- 
tor. The color image represents the out of plane component of 
the spin textures. All three currents are azimuthally symmet- 
ric (a) The spin currents have oscillations on the length scale 
r ~ l/v. On the axis with z — 0, the local spin vector is po- 
larized entirely out of plane at solutions to tanm = —rv/v 2 . 
The odd solutions feature in-plane vortex loops of spin, while 
the spin forms anti-vortices at the even solutions, (b) The 
orbitals are dominated by the in-plane component, which is 
stronger than the out of plane components by a factor of v. 
A small out of plane component is largest on the z = axis. 
In the upper half plane with y > 0, current converges on the 
point r — 0, while on the lower half plane all the current di- 
verges away from the point r = 0. (c) The total current is the 
sum of orbital and spin currents. The total, spin and orbital 
currents are all independently conserved. 



is the momentum-space current density. We find that the 
ground state current in momentum space is given by 



-(p-v) 2 



±- 



47T 



(31) 



for a state with m = ± h. 



2. Position Space 

In position space, the spin-orbit coupling is imaginary, 
so the continuity equation will have contributions from 
the orbital and spin degrees of freedom: 

j=Jo+j s (32) 

where the orbital current j = 3 WV^l is analogous 
to the orbital current in momentum space space, and 



7 



the spin current is j s = v [ip'aaj. For an eigenstate imuthally symmetric. The ground state currents can be 
of the trapped Weyl Hamiltonian, these currents are az- calculated exactly for large spin-orbit coupling, and are 

given by 



Jo = 
j s = 



1 

4tt 

v 

4tt 



(/o + 9r/«T - /<Wo + ) cos^f + sin^ (f-§ + f+$ 



(/+) 2 cos0-(/ o -r 



sin o r 



(fo + 



(fo 



cos I 



6+(2f+foSm6)4> 



(33) 
(34) 



for the ground state with j = i and m — \ . These cur- 
rents are shown in Figj4] For a plane defined by fixed 
<j), both currents are azimuthally invariant: that is the 
current fields do not change forms under a rotation by 
4>. Since the length scale \/v appears only in the radial 
variable, the radial component of the orbital current will 
be stronger than the polar and azimuthal component by 
a factor of v. As is seen in Fig. |4jb), the out of plane 
component of spin is negligible only near the plane de- 
fined by 9 = 7r/2. The lower half-plane has all currents 
flowing to the point r = 0, while all currents in the upper 
half plane diverge from this point. 

The spin currents do have an out of plane, or az- 
imuthal. component of spin. The in-plane component has 
a similar structure to the orbital current, in the lower 
half-plane the current converges to r — 0, while it di- 
verges from the point r = in the upper half plane. The 
out-of-plane component oscillates on the scale r ~ 1/v 
with amplitude cos 9. At special points corresponding to 
the solutions of tanru = —rv/v 2 , the out of plane spin 
current is completely polarized in the azimuthal direc- 
tion. The sum of the two currents is seen in Fig. |4jc). 
Since the orbital and spin currents are independently con- 
served, the total current is also conserved. 

IV. CONCLUSION 

The problem of the three-dimensional spin-orbit cou- 
pling in a harmonic trap was considered. The system 
has a conserved total angular momentum that allows the 
Hamiltonian to be projected into sectors where j, m are 
good quantum numbers. Each sector is tri-diagonal, pro- 
viding for an efficient numerical calculation of the spec- 
trum. At large values of the spin-orbit coupling parame- 
ter, the system undergoes a dimensional reduction from 
three to one. The spectrum is well approximated by 
Landau levels in the radial coordinate, with splittings in- 
versely proportional the the spin-orbit coupling strength. 

The conservation of total angular momentum allows 
us to express Schrodinger's equation as a set of coupled 
differential equations in the radial coordinate. We find 



the form of these differential equations in both position 
and momentum space. Using asymptotic analysis on the 
momentum space differential equations, we reproduce an 
analytic spectrum that well approximates the spectrum 
found numerically. We further use this approximation to 
find analytical expressions for the low energy eigcnfunc- 
tions in both position and momentum space. 

Finally, we explore the properties of the ground state 
wavefunction. We find a momentum space orbital cur- 
rent in the azimuthal direction. In position space the 
total current can be decomposed into a spin current and 
an orbital current. All currents are invariant under az- 
imuthal rotation. At large spin-orbit coupling, the spin 
currents have an azimuthal component that oscillates on 
the inverse of the spin-orbit strength. Along the z = 
axis, the spin currents alternate between complete polar- 
ization with and against the azimuthal plane. At these 
special points, the in plane spin density has vortex and 
anti-vortex structure. The orbital current is character- 
ized by a convergence of all currents at the point r = 
in the upper half-sphere, with a divergence of away from 
the same point in the lower half-sphere. 

A single trapped particle with Weyl coupling has rich 
ground state textures. This suggests that a system of 
trapped bosonic atoms with Weyl coupling will have 
novel many-body phases. Early evidence shows such a 
system will have qualitatively new ground states, such as 
a cubic lattice phase. |33ff55] 
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VI. APPENDIX A 



In this appendix we show that A + and A~ have the matrix elements as described by Eqj7]and|8] To do this, we first 
need to calculate the matrix elements of the three-dimensional isotropic harmonic oscillator basis, (n! V m! \a q \nlm) 
and (n' V m! \a^ q \nlm) , where q — 0, ±1, and the states \nlm) are the states three-dimensional harmonic oscillator basis 
as defined in the main text. 

It is convenient to define the number basis \n, I, mi) where (n, I, mi) are respectively the radial, angular and magnetic 
quantum numbers. In the absence of spin-orbit coupling we have the harmonic oscillator energies E = kuj (2n + l+ |) 
where a state of angular momentum I has a 21 + 1-fold degeneracy of mi . To do this we define the spherical creation 
and annihilation operators a± — {o-x T ia-y), a± = T^tj ( a x i ^ a y)> anc ^ a a — a z, o^o = o,\- We will find it helpful 
to express the angular momentum generators as 



L + 



\/2 {a) + a z + aja_j , 
V2 (a?_a z + a\a + ^j , 



a a_ 



(35) 
(36) 
(37) 
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with L x = \ (L + + L_) and L y = — | (L + — L_). If we additionally define the operators 

S+= 1 -(aiy -a{al (38) 
S_ = l - (a ) 2 - a+a_ (39) 

So = 1(^+1) (40) 

These operators commute with the angular momentum operators, and satisfy the commutation relations [S'+jS 1 -] = 
—2Sq and [So,5±] = ±£±- The three-dimensional harmonic oscillator eigenstates can be expressed as 

SlLt m \0,l,l) (41) 

where \0,l,l) = ( a +J |0, 0, 0) is a state of maximum angular momentum. It can be seen that a state \nlm) has 
energy E = 2n + I + §. 



\nlm) 



(l + m)\ 



i!r(n + i+f) (21)1(1 -m) 



A. Caculation of: a\ 

To calculate aL, we first note that the operator changes the total energy number by An = +1 unit, so it 
must connect states with either An = +1 and AZ = —1, or An = and Al = +1. Other states with An > 1 and 
Al = — 1 + 2(An — 1) are consistent with this condition, but must have zero matrix elements from the fact that a q 
is a spherical tensor of rank-1, and higher order corrections are inconsistent with angular momentum conservation. 
Using this fact, the matrix elements can be decomposed as 



(n', I', m'\ at \n, I, m) = (n, I + 1, m - 1| a[_ \n, I, m) £n',n^',/+i<W,m-i 

+ (n+ 1,1 - l,m - 1| a[_ \n,l,m) 5 n > .n+i^v ,i-i5 m > , m -i 



(42) 
(43) 



The two non-zero elements can be calculated individually. We first calculate (n, I + 1, m — 1| a)_ \n, I, m). It is conve- 
nient to use the equivalent definition of |n, I, to), given by 



\nlm) = 



r(/ + 1) (/-to)! 

iT {n + l+D (21)1(1 + m)\ 



SlL l + m \Q 1 l,-l) 



(44) 



to express the desired element as 



(n,Z + l,TO-l|a T _|n,Z,TO) = Af n , +hm ^Af n ,, m (0| (a_)' +1 (L_)'+ m (5_) n al (S + ) n (L + ) l+m (al) 1 |0) (45) 



= tf n ,i+i, m -itf n ,i, m (0| (a_y +1 (L_y+ m (S_nS + ral(L + ) l + m (aly \0) (46) 



where the normalization constant is Af n , 



r('+|) 



Q-mM 



and we have used the commutation rela- 



n\T(n+l+l ) (2l)\(l+ m )\ 

tion [ai,S + ] = 0. The bra (0| (a-) l+1 (L-) l+m is an eigenbra of the the operator (#_)" (S+) n with eigenvalue 
r^+3/2+i)" 1 ^ • ^ e can a S am use the commutation of the operator [L_,a[_] = to find the remaining factor 



(0| (a-) l+1 (L-y +m al(L + y +m (al) 1 |0) = (0| (o_) ,+1 ol (L_)' +m (L + )'+ m (a t _) ( |0) 

= (l + l) (0| (a-) l (L_y +m (L + ) l + m (aly |0) 

These factors combine to give the non-zero value of the matrix element 



(47) 
(48) 
(49) 



(n, I + 1, m — 1| a_ |n, l,m) = - 



1 /n + / + 3/2(Z-m + 2)(Z-TO+l) 



1 + 3/2 



1 + 1/2 



(50) 
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Similar techniques can applied to evaluate the other non-zero matrix element of the operator aL. The calculation is 
straightforward using the same techniques. The full matrix element is given by 



/ ' 'it, , \ 1 n + l + 3/2(l-m + 2)(l-m + l) 
(n ,1 ,m\a'_\n,l,m) = - J ^ + 3 / 2 1 + 1/2 °n> , n 6i' ,i+iO m > ,m-i (51) 



1 n + l (l + m)(l + m-l) 
~ 2]j 1 + 1/2 T^l/2 V,n+i*l',J-i<ym', m -i. (52) 

The remaining matrix elements can be calculated in an analogous way. The results can be summarized by the 
expression 

a q \n, I, m) = c q (n, I, m) \n, I + 1, m + q) + d q (n, I, m) \n + 1, 1 — 1, m + q) (53) 
where q = — 1, 0, 1, and the matrix elements are given by 



1_\ 1+191 / n + l + 1/2 
V2j \l (l + l/2)(l-l/2) 
i+kl 



c q (n,l,m) = ( ) \l n , , /oW> — r^\fq( l ' m ) ( 54 ) 



j \ i+H / n + Z + 3/2 
c+(n,Z,m) = ( -^J J {l + 3m + 1/2) 9*W ( 56 ) 



1 \ i+kl / 1 
1 \ 11 / n + l 



<(M,m) = (-!)■ ( -^j J _____ / _, (Z>m ) (57) 



where the functions f q (l,m) and g q (l,m) are defined as 



v / (/ + m)(Z + m- 1) g = +1 
f q (l,m)={ y/(l + m)(l-m) _q = (58) 
yj(l-m){l-m-l) q = -l 

y/(l + m + 2)(l + m + l) q = +l 
g q (l,m)=\ y/(l + m+l)il-m + l) q = ■ (59) 
y/(l-m + 2)(l-m + l) q = -l 



The matrix elements of the operators a q can be calculated through conjugation. 



B. Matrix elements of A + and A . 

We can now calculate the matrix elements of A + and A~ for states of good total angular momentum labeled by 
quantum numbers j, m. Recall in the main text that these states are defined by 




\n,j,m,\)=e t +* \ \ ' "\ (60) 

V V 20'+^) F'J+2' K 

where is an arbitrary phase. We will find it convenient to express this as 

r il hi + h m ~ 5) 

where 7 J = Ae*> ^/^^ and ^ = -e»> y/^^. 

We now consider the action of the operator er • p on the states \n,j, m, A), 



nJ,rn,X)=(i\ n >. + f m -\y), 



11 



(63) 



r 2 {(J _ a _ - a+a+) + « za \ \ n ,j, m , A) = {-^t\ IV m + ^ + r | |n '.^' m ) ^ 

V V27>- |n,j + f ,m- ±) - 7^ |n,j + f ,m + |) / 

(-\/2 7 j; c+(n, j + f ,m + i) +7lco (n, j + |,m - §)) |n, j - 1 + §,m - §) 

(\/27lc:(n,j + |,m- i) - 7 j;co (n,i+ |,m+ i)) |n,i - 1 + |,m+ i) 

(-^7j;d;(n,i + |,m+i)+7j;do(n,i + f,m- i)) |n - l,j + 1 + |,m - §\ 
(v^yldl(n, j + |,m- |) - 7t<^(n,i + im+ I)) |n - l,j + 1 + f,m + \) ' 

This can be expressed as 
where the coefficients 



n + j + | + 5 ( t ( ■ A 1 



2(i+| + |)(j+ 2 2 , 



/+ (j + 'i,m+ ^ ) 7i + /o (j + ^,m - ^ ) 7 J ) (67) 



" % + +xly{j + a^-\) ^ {VU + ^m + x x )(j + m + x x ) + - Am + x x ) (j - m + x x fy 



and 



ri(A)= (v^ci ^„,j + ^, m _0 _ 7 j;c- ^„,j + ^ ro + i^ (69) 



i0 A / (n + j + a; A )(j-m + a; A -l) 
6 V 4(. 7 +x A ) 2 (j+ a ; A -l) 

Explicitly calculating these for A = ±1 we get: 



j - Am + ai A ) (j - m + x A ) + \/ (j + Am + x x ) (j + m + x x Wl) 



V2~fic+ (n, j + \ , m + + 7W (n, j + \ , m - ) = ^2(n + j + l)^±^ (72) 
V2-yic+ ( n,j - *,m + i ) +7tc(7 ( n,j - *,m- * ) ] = (73) 



and 



V2 7 ; C I [nJ+^m--j - 7 ico (n, j + - , m + - j j = e*+ ^(n + j + l)^ (74) 
\/2^_cZ (n,j - \,m - \ ) - 7-^ (n,j - ^,m+ ^ ] ) = 0. (75) 



2 2 / V 2 2 

Similarly, we calculate the equivalent terms for T^(A) with s =t,|- This gives 



' n (j — m + x x + 1) 
Mj + x x ) 2 (j + x x + l) 



(76) 

v/ (j + Am + xx) (j - m + x x ) - \^ (j - Am + x x ) {j + m + x A )) (77) 
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2' m+ 2 



(78) 



= _ gt A , / n (j + m + za + 1) 

Explicitly, 



4(j + x A ^(j + XA + l) (aVC? - Am + x A ) (j +m + x x ) - y/(j - m + x x ) (j + Am + x A )) (79) 



rJ(+)= (-v^7t^ ("'•?' + ^> m + +^+ d o ( n 'J' + ^' m -0) =0 ( 8 °) 

^(-) = (nj-^m+l) +1 U- (nj-l,m-l)) =e »-^Ji^±± (81) 



ri(+)= [VZildZ [n,3 + \,m-^\ -7t4 + ) =° ( 82 ) 



ri(-) = (W<T («, j - - - 7 £< 



1 „ , !\\ ^-.rKZ lj + m + 1 



^-- r m + - 2 ))=e^-V2nJ^ TW (83) 



We summarize these relations as (restoring the subscript on n r ) 

A~ \n r ,j,m,+) = e lA ^^/2rh r \n r - l,j,m, -) (84) 
A' \n r ,j,m, -) = e- iA *V%+i + 1) \n r ,j,m, +) (85) 

where A<p = (<f> + — </>_). We can find the action of the operator A + = A~^ through conjugation 



A + \n r ,j,m,+) = e- tA<t > x /2(n r +j + 1) \n r ,j,m,-) (86) 
A + \n r ,j,m,-) = e lA < t > x /2(n r + 1) \n r + l,j, m, +) (87) 



We can show that \{A + ,A~) = N + 3/2, so 



H=±{A+,A-} + t ^=(A+-A-). (88) 
If we chose e jAl * = —1, we see that the matrix elements for p • a = (A + — A~) are: 

(ri r ,j,m,+\-^= (A + -A~) \n r ,j,m,-) = i (\Av + lS nr+liK - y/n r +j + l<5„ r ,„^ (89) 

(ri r ,j,m,-\-j= (A + -A~) \n r ,j,m,+) = i (y/n r + j + 16 nrtn > r - V^V-ix) ■ ( 90 ) 

VII. APPENDIX B 

In this appendix we calculate the radial Schrodinger equation in position space by Fourier transforming the 
momentum-space version. We want to find the Fourier transform of the eigenfunctions of the harmonic oscilla- 
tor. We first need to find the momentum space eigenfunctions of the three-dimensional spherical harmonic oscillator 
int terms of the momentum space harmonic oscillator wavefunctions(p|^„; TO ) = f , dr e lp ' r (r\ip n i m ) . Using the 

V(2f) 3 

expansion of plane waves into spherical harmonics, e tp r = 47r^ im i l ji{pr)Y l m (f)) (F; m (f))*, we find the relation 

(p|Vw) = / 7 |==4 7 r^ 'jvlprW '(p) (V (r))* (r#„ (m ) (91) 

VI / l' m ' 

^YT(P) I dr {r 2 3l (pr)R ln {r)) ■ (92) 
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To find the radial Schrodinger's equation for the Weyl coupling in position space, 

1pnjm{p) = E nj 1p njm (p) 



V 2 , P 2 



2 +Y + vp - (T 



(93) 



it is easiest to begin with the same equation in momentum space, and Fourier transform the eigenstates ip n j m (p) = 
J e tp r tfi n j m (r)dr. But recall, we have showed that j and m are good quantum numbers, and the eigenfunctions have 

the form Vv, m (p) = fn(p)xt m (p) + fn(p)xj m (p)- Tne spinors xf m contain spherical harmonics of order I = j =F \, 
so they are eigenfunctions of the Fourier transform. The position space wavefunction thus has the form 



<f>njm{r) = 4tt y 



d 3 r 



J2(-i?Mpr)Yr(p) on*))* 



_ lm 



{fn(p)Xj m (p) + fn (p)X jm (P)) 



(94) 
(95) 



where g±(r) = J p 2 ji (pr) f±(p)dp, and I = j =p \ - 

We are now ready to transform the momentum-space Schrodinger equation, we first multiply by e 4p r , and then 
integrate over momentum to get 



rf3 P e i P - 



K , p 2 



VVjm(p) = 



nj 



I 



d 3 p 



=e* p - r r/'„jm(p)- 



Transforming this term-by-term, the momentum space kinetic term becomes a position space trapping term 



/ 



Similarly, the momentum space trap becomes a kinetic energy term in position space 



v 2 



The spin-orbit term is more complicated, we first use the property that p • <rxf m (p) = PxJ m (p) 

I [e ip - r v(p-cr)^ njm (p)] =v J -J^L= [e < P-'p(/+(p)x7 rn (p) + / n -(p)x+,(p))] 

we then expand the exponential e lp r in terms of spherical harmonics, 



(96) 

(97) 

(98) 
(99) 

(100) 



d 3 i 



[e^v(p-a)^ njm (p)]=v-^=^i l J dp J dpp 3 [j t (pr)lT (?) (*T (p))* (/+(p)x7 to (p) + /» (p)x+ m ((p-P>]) 

* / dp(pVi(P r )/n J -(p))x+,(?) + / *(f 3 i J -i(pO/+(p))x7 m O:2) 



The spinors Xj" m have angular momentum components 2 = j =p | . These angular momentum variables do not match 
the angular momentum corresponding to the spinors multiplying them. This, along with the extra factor of p in 
the integrand J dp (p 3 ji(pr)f(p)) = J dpp 2 (pji(pr)f(p)), suggest it would be helpful to express the term pjj 1 (pr) in 

terms of spherical Bessel functions of the order We use the identities 



zjl(z) = z ^3l+i(z) + ( l + 2 )ii+i( z ) 
zji(z) = + ( l - l )ji-i( z ) 



(103) 
(104) 
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to express 



pji{pr) = 



d 1 + 2 
dr r 



■ _A l ~ 1 



3i+i(pr) 
ji-i(pr). 



(105) 
(106) 



These relations convert the product pji(pr) into a differential operator in r that can be removed from the integral. 
Together, these give the expression for the Fourier transform of the three-dimensional spin-orbit coupling term 



/ 



d 3 p 



[e ip ' r v(p ■ a) Vv?m(p)] = v 



+ 



TT 

dr r 



d J- 2 
dr r 



J dp (p 2 j._i (pr)/".(p)) x+ m (r) (107) 



(108) 



= w- 7 2 



+ — 

dr r 



3 + 



(109) 



Finally, we combine these to get the expression for the radial Schrodinger equation for the three-dimensional spin-orbit 
coupling in position space, 



1 / d 2 



dr 2 
1 / d 2 



r + 



r \ dr 2 



(j+M+j) r2 



r| + (i±M^i) +r2 



, d .7 + 9 



_d_ 
dr 



9ni( r ) = E nj9nj( r 



9ni( r ) = E nj9 nj (r 



(110) 
(111) 



